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ABSTRACT 

We  provide  new  results  enabling  robust  interferometric  image  reconstruction  in  the 
presence  of  unknown  aperture  piston  variation  via  the  technique  of  Redundant  Spac¬ 
ing  Calibration  (RSC).  The  RSC  technique  uses  redundant  measurements  of  the  same 
interferometric  baseline  with  different  pairs  of  apertures  to  reveal  the  piston  variation 
among  these  pairs.  In  both  optical  and  radio  interferometry,  the  presence  of  phase¬ 
wrapping  ambiguities  in  the  measurements  is  a  fundamental  issue  that  needs  to  be 
addressed  for  reliable  image  reconstruction.  In  this  paper,  we  show  that  these  ambigu¬ 
ities  affect  recently-developed  RSC  phasor-based  reconstruction  approaches  operating 
on  the  complex  visibilities,  as  well  as  traditional  phase-based  approaches  operating  on 
their  logarithm.  We  also  derive  new  sufficient  conditions  for  an  interferometric  array 
to  be  immune  to  these  ambiguities  in  the  sense  that  their  effect  can  be  rendered  benign 
in  image  reconstruction.  This  property,  which  we  call  wrap-invariance ,  has  implica¬ 
tions  for  the  reliability  of  imaging  via  classical  three-baseline  phase  closures  as  well 
as  generalized  closures.  We  show  that  wrap- invariance  is  conferred  upon  arrays  whose 
interferometric  graph  satisfies  a  certain  cycle-free  condition.  For  cases  in  which  this 
condition  is  not  satisfied,  a  simple  algorithm  is  provided  for  identifying  those  graph 
cycles  which  prevent  its  satisfaction.  We  apply  this  algorithm  to  diagnose  and  correct 
a  member  of  a  pattern  family  popular  in  the  literature. 

Key  words:  techniques:  interferometric  -  techniques:  image  processing  -  atmospheric 
effects  -  methods:  analytical  -  methods:  numerical 


1  INTRODUCTION 

Optical  interferometry  is  a  multi-aperture  imaging  technique 
which  is  attracting  increasing  interest  in  the  astronomical 
and  remote-sensing  communities.  The  appeal  of  this  tech¬ 
nique  is  primarily  due  to  the  high  resolution  it  affords  rela¬ 
tive  to  single- aperture  imaging.  Namely,  the  angular  resolu¬ 
tion  of  a  single  aperture  is  limited  by  diffraction  to  jy ,  where 
A  is  the  wavelength  of  the  light,  and  D  is  the  diameter  of 
the  aperture.  On  the  other  hand,  the  achievable  angular  res¬ 
olution  of  an  interferometer  is  instead  given  by  g-d — ,  where 
Bmax  is  the  maximum  spatial  separation  of  any  two  tele¬ 
scopes  in  the  array.  Therefore  with  interferometry  one  can 
achieve  the  same  high  resolution  offered  by  an  extremely 
large  (and  often  prohibitively-costly)  telescope  by  interfer¬ 
ing  light  from  several  telescopes  of  practical  size.  Optical 
interferometers  image  a  scene  by  sampling  the  2D  Fourier 
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Transform  of  the  scene.  Several  excellent  surveys  exist  which 
describe  the  concept  of  interferometry,  including  Labeyrie 
et  al.  (2006),  Glindemann  (2011).  Each  pair  of  telescopes 
measures  a  single  angular  spatial  frequency  of  radians, 
where  b  is  the  vector  difference  of  the  telescope  positions, 
which  is  known  as  a  baseline.  For  an  array  of  N  apertures, 
the  data  set  then  consists  of  all  such  measurements. 

A  principal  challenge  in  interferometry  is  variation  in 
the  complex  gains  among  the  multiple  apertures  of  the  in¬ 
terferometer.  In  radio  interferometry,  this  variation  can  arise 
from  differences  among  the  analog  components  of  the  an¬ 
tenna  elements  (e.g.  cable  length  differences)  in  the  array. 
In  optical  interferometry,  atmospheric  turbulence  distorts 
the  wavefronts  arriving  at  each  telescope  aperture  so  that 
their  effective  path  lengths  from  the  target  (or  optical  pis¬ 
tons)  are  altered  by  a  random,  non-uniform  amount.  For 
simplicity,  we  will  heretofore  refer  to  such  aperture-specific 
phase  variation  as  optical  path  difference  (OPD),  regardless 
of  its  source.  As  a  result  of  OPD,  the  Fourier  component 
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measured  by  apertures  i  and  j  is  given  by  yij  =  \yij\ej0iJ , 
in  which  §ij  =  (6ij  +  Sfiij)  mod  2n  1 ,  and  Gy  is  the  true 

Fourier  phase  at  spatial  frequency  and  dfiij  represents 
the  OPD  observed  between  apertures  i  and  j.  One  approach 
to  eliminate  OPD  is  to  form  triple  products  of  the  Fourier 
components  along  the  sides  of  a  baseline  triangle  (e.g.  b  12, 
623,  and  631).  Note  that  the  OPD  cancels  in  these  products 
and  hence,  like  the  Fourier  magnitudes,  these  so-called  bis¬ 
pectra  are  OPD-invariant  observables.  However,  for  a  non- 
redundant  array  with  distinct  baselines,  recovery  of  the 
Fourier  phases  from  the  bispectra  phases  is  an  ill-posed  prob¬ 
lem  since  there  are  only  (A^~1)  independent  bispectra  (Read- 
head  et  al.  1988).  Successful  bispectra-based  image  recon¬ 
struction  remains  feasible  in  spite  of  this  ill-posedness  (see 
e.g.  Thiebaut  (2013),  Besnerais  et  al.  (2008)),  but  in  this 
case  prior  constraints  (e.g.  on  the  image  support)  must  be 
enforced  to  regularize  the  reconstruction. 

An  alternative  and  intrinsically  well-posed  approach 
to  prior-regularized  reconstruction  is  to  use  baseline  re¬ 
dundancy  to  explicitly  solve  for  OPD  variation;  an  array 
with  baseline  redundancy  contains  repeated  instances  of 
the  same  vector  baseline  involving  distinct  aperture  pairs. 
Since  Fourier  phases  can  be  assumed  to  be  equal  for  all 
repeated  vector  baselines,  an  observed  difference  amongst 
their  corresponding  measurements  exposes  the  contribution 
of  the  OPD.  This  idea  of  using  redundant  arrays  to  calibrate 
out  OPD  variation,  known  as  redundant  spacing  calibration 
( RSC ),  was  developed  in  works  such  as  those  by  Arnot  et  al. 
(1985)  and  Greenaway  (1990).  In  recent  years,  innovation 
in  optical  technology  has  engendered  a  revival  of  interest 
in  the  RSC  technique.  The  simultaneous  (or  Fizeau- style) 
measurement  of  fringes  on  a  common  focal  plane  has  long 
been  a  popular  method  of  acquiring  many  baseline  mea¬ 
surements  in  an  economical  manner.  However,  the  Fizeau 
method  had  been  incompatible  with  RSC  techniques  since 
the  fringes  formed  by  each  set  redundant  baselines  would 
alias  on  the  focal  plane.  An  elegant  solution  to  this  problem 
was  proposed  by  Perrin  et  al.  (2006).  This  work  developed 
the  idea  of  segmenting  the  entrance  pupil  of  a  single  tele¬ 
scope  into  an  RSC  arrangement  of  sub-pupils  from  which 
the  light  was  then  coupled  via  single  mode  fiber  to  a  non- 
redundant  exit  pupil,  thereby  permitting  unambiguous  and 
simultaneous  fringe  detection  for  an  RSC  array.  A  recon¬ 
struction  algorithm  for  this  architecture  was  then  proposed 
in  Lacour  et  al.  (2007).  Even  more  recently,  RSC  has  been 
implemented  as  the  calibration  scheme  of  choice  for  several 
radio  interferometers:  the  Donald  C.  Backer  Precision  Array 
for  Probing  the  Epoch  of  Reionization  (PAPER)  in  South 
Africa  (see  Ali  et  al.  (2015)),  the  MIT  Epoch  of  Reionization 
(MITeOR)  in  the  United  States  (see  Zheng  et  al.  (2014)), 
and  the  Ooty  Radio  Telescope  (ORT)  in  India  (see  Marthi 
&,  Chengalur  (2014)). 

As  will  be  shown  below,  N- 3  independent  redundant  re¬ 
lations  are  required  for  unique  determination  of  atmosphere 
and  Fourier  phases  -  a  condition  we  will  refer  to  as  critical 
redundancy.  An  oft-cited  drawback  of  the  RSC  approach  is 
that  it  reduces  the  number  of  unique  spatial  frequencies  mea¬ 
sured  by  the  interferometer.  However,  as  Figure  1  illustrates, 


1  In  this  paper,  we  use  mod  to  denote  the  modulo  operation 


Number  of  apertures  in  array 

Figure  1.  Fraction  of  redundant  baselines  required  for  critical 
redundancy  vs.  aperture  count 

the  fraction  of  distinct  uv-samples  sacrificed  for  critical  re¬ 
dundancy  becomes  increasingly  negligible  as  the  number  of 
apertures  in  the  array  increases.  Nevertheless  the  RSC  tech¬ 
nique  presents  other  challenges  which  must  be  overcome  for 
reliable  imaging.  Central  among  these  challenges  is  the  prob¬ 
lem  of  integer  phase  ambiguities  which  arise  from  the  fact 
that  the  interferometric  phase  is  only  known  modulo  2 n. 
Indeed  these  ambiguities  have  been  shown  to  play  an  im¬ 
portant  role  in  accurately  recovering  sensor  complex  gains 
and  object  complex  visibilities  while  imaging  with  real  inter¬ 
ferometric  instruments  (see  e.g.  Liu  et  al.  (2014),  Eastwood 
et  al.  (2009)).  In  this  paper,  we  describe  these  ambiguities 
and  how  they  can  be  mitigated  using  a  combination  of  lat¬ 
tice  theory  algorithms  and  careful  array  design.  We  will  see 
that  these  ambiguities  have  a  fundamental  presence;  namely, 
they  exist  whether  the  calibration  strategy  works  with  com¬ 
plex  visibilities  (which  we  call  the  Phasor  approach)  directly, 
or  their  respective  logarithms  (which  we  call  the  Phase  ap¬ 
proach).  To  the  best  of  our  knowledge,  the  results  in  this 
paper  are  the  first  to  provide  array  conditions  allowing  un¬ 
ambiguous  interferometric  phase  determination  in  spite  of 
wrap  ambiguities  in  the  Phase  approach,  and  corresponding 
false  minima  in  the  objective  of  the  Phasor  approach. 

To  motivate  the  analysis  in  the  paper,  we  provide  an  ex¬ 
ample  illustrating  the  effect  that  wrap  ambiguities  can  have 
in  RSC-based  image  reconstruction.  Consider  the  pattern 
depicted  in  the  Figure  2.  This  pattern  belongs  to  one  of  the 
more  popular  array  classes  in  the  interferometry  literature: 
the  so-called  Y-patterns  (see  e.g.  Arnot  et  al.  (1985),  Blan¬ 
chard  et  al.  (1996),  Labeyrie  et  al.  (2006),  Eastwood  et  al. 
(2009),  Liu  et  al.  (2014)).  The  corresponding  spatial,  or  UV , 
sampling  is  provided  in  the  right  panel  of  the  Figure. 

To  demonstrate  the  potential  effect  of  wrap  ambiguities 
on  reconstruction,  we  simulated  both  Phase-  and  Phasor- 
based  reconstruction  from  noiseless  measurements  with  this 
pattern.  The  results  are  shown  in  Figure  3.  The  upper  left 
panel  shows  the  true  image,  and  the  upper  right  panel  shows 
the  inverse  Fourier  transform  of  the  UV-sampled  visibility 
function  of  the  object  (i.e.  the  so-called  dirty,  or  interfero¬ 
metric  image).  The  lower  left  panel  shows  the  reconstruction 
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Figure  2.  Y-pattern  Array  Example 


Figure  3.  Reconstruction  Results  for  Y-Pattern  Example 


result  with  an  implementation  of  the  Phase  method  similar 
to  that  developed  by  Lannes  (2003),  and  the  lower  right 
panel  shows  the  same  for  an  implementation  of  the  Phasor 
method  (Marthi  V  Chengalur  2014)  (Lacour  et  al.  2007). 
Reconstruction  suffers  from  phase  wrapping  error  in  the  for¬ 
mer  case,  and  a  corresponding  false-minimum  trap  in  the 
Phasor  case.  In  the  course  of  this  paper,  we  will  first  iden¬ 
tify  this  ambiguity  from  a  mathematical  perspective,  relate 
it  to  a  particular  physical  structure  (i.e.  the  existence  of  a 
certain  type  of  loop  in  the  interferometric  graph) ,  and  pro¬ 
vide  a  simple  algorithm  for  identifying  such  structures  in  an 
arbitrary  array  so  that  they  can  be  remedied. 

The  paper  is  organized  as  follows.  In  Section  II,  we  re¬ 
view  previous  work  on  the  integer  ambiguity  problem,  and 
discuss  its  presence  in  the  Phase  approach.  We  provide  new 
mathematical  conditions  for  an  aperture  pattern  to  be  wrap- 
invariant ,  meaning  that  the  effect  of  the  2^-periodicity  of  the 
measured  interferometric  phase  can  be  eliminated  in  image 
reconstruction.  These  results  are  founded  upon  techniques 
from  lattice  theory,  as  well  as  the  well-known  Smith  Normal 
Form  (SNF)  of  an  integer  matrix.  We  show  the  implications 
of  these  results  on  imaging  with  three  types  of  interfero¬ 
metric  observables:  the  baseline  phase  measurements,  their 
traditional  closure  phases,  and  generalized  closure  phases.  In 
Section  III,  we  relate  these  mathematical  conditions  to  con¬ 
ditions  on  the  aperture  pattern  itself.  Namely  we  show  that 
wrap- invariance  is  conferred  upon  arrays  satisfying  a  certain 


loop-free  condition.  As  an  illustrative  example,  we  diagnose 
a  pattern  belonging  to  the  popular  Y -pattern  class  and  rem¬ 
edy  it  to  be  loop-free.  Finally,  in  Section  IV,  we  show  that 
the  computationally-complex  SNF-based  approach  for  ambi¬ 
guity  resolution  is  actually  not  necessary  for  wrap- invariance 
in  many  cases,  as  long  as  a  wrap-induced  image  shift  can  be 
tolerated  2.  For  such  scenarios,  we  show  that  wrap- invariance 
provides  a  certificate  for  the  success  of  existing  Phase  and 
Phasor-based  approaches  in  avoiding  wrap-induced  errors  in 
the  former,  and  false  global  minima  in  the  latter.  Finally  we 
summarize  our  results  in  Section  V. 


2  PHASE  WRAPPING  AMBIGUITIES  IN  RSC 
IMAGE  RECONSTRUCTION 

In  this  Section,  we  describe  a  Phase  Approach  algorithm 
leveraging  the  CVP-approach  for  phase-wrap  resolution, 
which  to  the  best  of  our  knowledge  was  first  developed  in 
Lannes  &  Anterrieu  (1999).  We  begin  by  identifying  the 
fundamental  phase  ambiguity,  and  subsequently  assess  its 
impact  on  the  phase  error  in  the  RSC  phase  solution.  In 
the  process,  we  develop  mathematical  conditions  for  wrap- 
invariance  which  will  form  the  basis  for  the  notion  of  a  wrap- 
invariant  pattern  in  Section  3.  Finally,  we  assess  our  results 
in  the  context  of  similar  results  for  approaches  requiring  the 
use  of  closure  phases  (Lannes  2003). 


2.1  Identifying  the  Fundamental  Phase  Ambiguity 

The  traditional  approach  to  RSC  reconstruction  operates  on 
the  measured  baseline  phases  (see  e.g.  Arnot  et  al.  (1985), 
Greenaway  (1994)).  To  illustrate  the  approach,  let  us  con¬ 
sider  an  interferometer  which  operates  at  a  wavelength  A 
with  two  apertures  (say,  i  and  j)  separated  by  a  vector  base¬ 
line  distance  of  by.  In  the  absence  of  any  optical  path  differ¬ 
ence,  the  interference  pattern  formed  by  these  two  apertures 
encodes  a  sample  of  the  object’s  Fourier  Transform  at  spatial 
frequency  .  Let  the  true  Fourier  phase  (which  we  will  re¬ 
fer  to  as  object  phase),  measured  by  this  interference  pattern 
be  denoted  as  0[j.  The  measured  phase  is  given  by: 


Pij  —  @ij  +  (pj  (pi  +  Ane  (1) 

where  (pj-(pi  is  the  optical  path  difference  between  apertures 
j  and  /,  and  e  is  unknown  phase  wrap  integer  arising  from 
the  fact  that  interferometric  phase  measurements  are  only 
known  modulo  2/r. 

Consider  an  interferometric  array  which  simultaneously 
makes  many  such  measurements  amongst  its  N  apertures. 
Suppose  that  of  the  array’s  baselines,  d  of  them  are 
distinct.  Further  suppose  we  have  a  solution  set  {(pi)  and 
{Oij}  for  these  equations.  Let  ry  denote  the  vector  position 
of  the  i-th  aperture.  As  noted  by  several  authors  (see,  e.g. 
Wieringa  (1992)),  we  can  obtain  another  valid  solution  set 
simply  by  replacing  each  0*  with  0^  =  0;  +  0q  +  z  •  r*,  and 
each  Oij  with  6?.  =  6 ij  -  z  •  (rj  -  Vf),  for  arbitrary  0q  and 


2  As  will  be  shown,  this  image  shift  is  distinct  from  the  fundamen¬ 
tal  degeneracy  of  object  position  in  interferometric  measurements 
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Figure  4.  Example:  5-aperture  RSC  pattern.  The  six  distinct 
vector  baselines  are  shown. 


z.  Since  the  free  vector  z  is  a  two-parameter  vector  rep¬ 
resenting  the  inherently-ambiguous  position  of  the  image 
within  the  Field-of-View  and  the  free  parameter  0q  is  sim- 
ply  a  scalar  piston  offset,  the  kernel  of  the  RSC  system  is 
three-dimensional.  This  is  the  tilt-position  degeneracy  which 
is  fundamental  in  interferometric  reconstruction.  The  end 
result  is  that  the  RSC  system  contains  d  unknown  distinct 
object  phases  and  N  unknown  aperture  pistons,  and  is  rank- 
deficient  by  at  least  3.  This  implies  that  there  are  at  most 
(d  +  N  -  3)  linearly-independent  equations  in  the  RSC  sys¬ 
tem,  and  hence  at  most  N  -  3  redundant  relations  can  be 
linearly-independent.  We  will  assume  for  the  remainder  of 
the  paper  that  our  array  contains  N  -  3  independent  rela¬ 
tions.  Under  this  assumption,  we  could  in  principle  solve  for 
a  particular  solution  of  this  system  by  arbitrarily  setting  two 
object  phases  (whose  spatial  frequencies  are  not  co-linear) 
and  one  piston  phase.  This  particular  solution  would  then 
differ  from  the  true  solution  by  a  phase  ramp  in  the  Fourier 
domain,  corresponding  to  an  image  shift  in  the  spatial  do¬ 
main.  In  this  section  we  will  instead  construct  a  different 
particular  solution  which  is  immune  to  the  effects  of  phase¬ 
wrapping  by  design. 

As  an  example,  consider  the  simple  array  in  Figure  4. 
There  are  (2)  -  10  baselines,  of  which  4  are  redundant. 
A  critical  array  of  5  apertures  would  have  2  redundancies. 
Therefore  this  array  possesses  more  redundancies  than  nec¬ 
essary  (call  it  strongly  redundant ),  and  we  anticipate  that 
the  resulting  system  will  be  overdeter  mined. 

The  measurement  equations  associated  with  this  array 
can  be  written  in  matrix  form: 


1  0  0  0  0 

0  10  0  0 

0  0  10  0 

0  0  0  1  0 

0  0  0  0  1 

0  0  0  1  0 

0  0-100 
-1  0  0  0  0 

0  0  0  0  0 

0  10  0  0 


0  1-10 

0  0  1-1 

0  0  0  1 

0  0  0  0 

0  10-1 

0  0  10 

0  10  0 

0  0  0  1 

10  10 
0  10  0 


Denoting  the  matrix  above  as  M, 
systems  in  compact  form  as: 


M 


6 

0 


P  +  2/re 


0  0 

0  0 

-1  0 

1  -1 

0  0 

-1  0 

-1  0 

0  -1 
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045 
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025 

01 

02 

03 

04 

05 


=  P  +  2/re  (2) 


we  can  write  such 


(3) 


The  example  above  illustrates  that  the  general  phase 
measurement  matrix  will  have  two  sets  of  columns:  one  cor¬ 
responding  to  the  object  phases,  and  one  corresponding  to 
the  path  differences.  Adopting  the  notation  of  Lannes  &  An- 
terrieu  (1999),  let  K  denote  the  subspace  spanned  by  the  first 
set  of  columns,  and  L  the  subspace  spanned  by  the  second 
set. 

If  we  let  n  -  be  the  number  of  baselines  in  the 
array,  the  phase  measurement  matrix  M  will  be  of  size  n  X 
(d  +  N) .  For  a  strongly-redundant  array  like  the  one  in  the 
example  above,  the  column-space  K  +  L  of  the  matrix  will 
clearly  not  span  Rn.  Therefore  the  wrapped  measurement 
vector  p  will  not  in  general  fall  in  the  the  subspace  K  +  L 
(and  potentially  can  be  quite  far  from  it).  In  the  absence 
of  measurement  noise,  we  can  unwrap  these  measurements 
by  identifying  those  integer  correction  vectors  e  for  which 
P*  =  P  +  2ne  lies  in  K  +  L.  In  the  presence  of  noise,  on  the 
other  hand,  the  unwrapped  vector  will  not  generally  lie  in 
K  +  L  (but  for  low-to-moderate  noise  will  be  in  the  vicinity). 
Thus  we  search  for  vector(s)  fi*  which  are  as  close  to  K  + 
L  as  possible  in  a  weighted  least-squares  sense  (Lannes  & 
Anterrieu  1999),  i.e.  we  search  for  the  vector 


w(/T(e)-M 

(4) 

where  W  is  the  weighting  matrix.  If  we  let  L  denote  the 
phase  measurement  covariance  matrix  and  set  W  =  L_1,  this 
is  equivalent  to  searching  for  vectors  e  which  minimize  the 
projection  of  a  whitened  measurement  W2  p*  =  W2  (p+2ne) 

onto  the  space  (K  +  L)^  :=  ker((\N 2M)  ),  where  ker( A) 
denotes  the  kernel,  or  nullspace,  of  the  matix  A.  Specifically 
we  seek  to  minimize: 

/(e)  =  ||PwW2  (/5  +  2/re)||2  (5) 

where  pW  is  a  matrix  representing  the  orthogonal  projection 
from  Rn  onto  ( K  +  L)^. 

Letting  e'  =  -e,  we  can  rewrite  the  above  objective 
function  as: 

/(e')  =  ||PwW2(/?-2/re')l|2  =  ||PwW2/5  -  2/rPwWi  e’\\ 

(6) 

Lannes  &  Anterrieu  (1999)  showed  that  this  optimiza¬ 
tion  problem  is  equivalent  to  the  so-called  closest  vector 
problem  in  the  theory  of  lattices.  We  will  define  a  lattice 
L(Zn)  as  the  set  of  points  generated  by  integer  combinations 
of  the  column  vectors  of  a  matrix  L.  Letting  P  =  P^W2,  our 
optimization  problem  then  amounts  to  the  following:  Find 
the  lattice  point  in  P(Zn)  which  is  closest  to  Pp/(2n).  A 
compact  representation  of  the  lattice  T  is  given  by: 


6rsc 

rRSC  = 

tpRSC 

=  argmine^,(f, 

Im<n-(d+N-3) 

L  aivi  I  Va«  e  : 

i= 1 


(7) 


where  {i7j}  are  linearly-independent  and  together  form  a  ba¬ 
sis  of  the  lattice.  Given  the  lattice  basis,  several  algorithms 
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exist  for  finding  the  closest  lattice  point  to  a  specified  vec¬ 
tor.  A  popular  class  of  algorithms,  known  as  the  Sphere- 
Decoding  algorithms,  are  efficient  searches  for  the  closest  lat¬ 
tice  point  within  a  hypersphere  of  a  certain  radius  centered 
on  the  input  vector  (see  e.g.  Agrell  et  al.  (2002)).  For  the  sim¬ 
ulations  in  this  paper,  we  instead  use  the  lower-complexity 
Babai  Nearest  Plane  (Babai-NP)  algorithm  (Babai  1986). 
For  lattice  bases  which  are  nearly  orthogonal  (such  as  those 
we  use  for  our  simulations),  this  algorithm  offers  reliable, 
albeit  not  guaranteed,  performance  in  practice. 

Suppose  we  have  found  a  basis  for  the  lattice  P(Zn), 
and  we  have  solved  the  Closest  Vector  Problem  for  a  given 
measurement  vector  yd.  Let  b*  be  the  output  of  the  Babai 
Nearest  Plane  Algorithm  -  i.e.  it  is  the  lattice  point  which 
is  the  closest  to  p.  We  now  seek  to  solve  for  the  wrap  vector 
corresponding  to  this  lattice  point,  i.e.  we  seek  a  solution  to: 

6*  =  Pe  (8) 

Note  that  P  is  a  (weighted)  projection  matrix  and  thus 
not  full-rank,  and  therefore  there  will  be  infmitely-many  so¬ 
lutions  to  this  equation.  The  Closest- Vector-Problem  algo¬ 
rithm  will  provide  one  particular  solution  ep .  The  complete 
set  of  solutions  is  then  given  by: 


K  +  L? .  In  either  case,  the  choice  of  the  residual-minimizing 
vector  yd*  is  not  unique.  To  see  this,  let  the  smallest  possible 
residual  norm  among  all  unwrapped  candidates  be  denoted 
as  rmin  •  The  set  of  unwrapped  measurement  vectors  rmin 
away  from  K  +  L  can  be  represented  as  discrete  points  in  a 
plane  parallel  to  K  +  L,  each  of  which  corresponds  to  a  dis¬ 
tinct  choice  of  ep .  Within  this  family,  consider  the  optimum 
vector  yd*  whose  corresponding  least-squares  solution  trsc 
has  minimal  error.  Any  choice  of  an  unwrapped  vector  in 
Equation  (10)  is  within  an  error  vector  2/re^  of  yd^,  where 
is  an  integer  vector  in  K  +  L.  RSC  algorithms  are  fun¬ 
damentally  blind  to  such  errors;  distinct  unwrapped  vectors 
yd*  and  ydg  both  produce  solutions  to  Equation  (4)  in  the 
noiseless  case,  as  do  their  respective  projections  onto  K  +  L, 
PK+L  and  Pqk+L'  noisy  case.  The  property  of  wrap- 

invariance  introduced  in  this  paper  ensures  that  the  effect 
of  such  an  error  on  the  phase  of  the  resulting  RSC  solution 
trsc  is  either:  merely  a  benign  integer  multiple  of  2 n  (c.f. 
Section  2.2),  or  a  linear  gradient  in  the  estimated  Fourier 
phases,  which  is  equivalent  up  to  an  extra  image  shift  (c.f. 
Section  4)  in  the  reconstructed  image.  In  order  to  develop 
conditions  for  wrap- invariance  as  they  relate  to  pattern  de¬ 
sign,  we  must  first  characterize  the  effect  of  the  residual  wrap 
error  on  the  RSC  least-squares  solution.  This  is  the  aim  of 
the  next  section. 


e  =  ep  +  eh  (9) 

where  e ^  is  any  integer  vector  in  the  kernel  of  P.  Suppose  we 
choose  one  such  vector  and  correct  our  phase  measure¬ 
ment  vector  accordingly.  The  corrected  phase  measurement 
vector  can  be  written  as: 

P*  =  P  +  2n(ep  +  cfc)  (10) 

Lemma  2.1:  e  K  +  L,  Ve^ 

Proof:  The  fact  that  e ^  e  &<?r(PwW2)  implies  that  W2  e ^  e 
ker(P\i\i).  This  in  turn  implies  that  W^e/j  e  im(\N^M)  and 
hence  that  e  im( M)  since  W2  is  invertible  by  construc¬ 
tion.  □ 

While  any  choice  of  e  K  +  L  will  admit  a  solution 
to  Equation  (4),  let  us  consider  the  optimal  one  e^o  which 
minimizes  the  error  in  the  ultimate  phase  solution  trsc- 
Let  ydg  =  yd  +  2n(ep  +  e^o)  be  the  corresponding,  optimal 
unwrapped  measurement  vector.  From  Lemma  2.1,  we  see 
that  the  unwrapped  vector  yd*  differs  by  some  2 ;re*  in  K  +  L 
from  ydQ,  i.e. 

0*  =  /3*o+  2ne*h  (11) 

The  impact  of  this  fundamental  ambiguity  is  the  main 
subject  of  this  paper.  We  have  depicted  the  situation  in  Fig¬ 
ure  5  to  provide  a  visual  summary  of  the  linear  algebra  in¬ 
volved.  As  we  have  seen,  the  RSC  process  begins  with  the  in¬ 
terferometric  phase  measurement  yd,  which  due  to  wrapping 
will  in  general  lie  far  from  the  range  K+L  of  the  measurement 
matrix.  By  solving  the  Closest- Vector-Problem  using  a  lat¬ 
tice  algorithm  such  as  the  Babai  Nearest  Plane  algorithm,  it 
is  possible  to  find  a  particular  correction  vector  2nep  which 
when  added  to  yd  minimizes  the  residual  in  Equation  (4), 
i.e.  the  weighted  distance  from  K  +  L.  In  the  noiseless  case, 
this  unwrapped  vector  yd*  will  lie  in  K  +  L  (i.e.  zero  residual), 
whereas  in  the  noisy  case,  it  will  in  general  remain  outside  of 


2.2  Quantifying  the  Effect  of  the  Fundamental 
Ambiguity 

Let  us  now  quantify  the  effect  of  this  unresolved  wrap  er¬ 
ror  on  the  ultimate  least-squares  solution,  which  can  be 
accomplished  in  two  easy  steps.  Following  standard  least- 
squares  principles,  we  first  find  the  projection  pK+L  of  the 
unwrapped  yd*  onto  K  +  L  whose  weighted  distance  from  yd* 
is  minimized.  Note  that  the  term  in  yd*  is  already  in  K  +  L 
and  hence  unchanged  by  projection.  Hence  we  obtain: 

0K+L  =  W4(l  -  PW)WM*  =  P*0,K+L  +  2neh  (12) 

where  Pq  K+L  is  the  projection  of  yd*  onto  K  +  L.  We  then 
solve  the  system: 


M 


;  Pk+l 


(13) 


As  aforementioned,  M  is  rank-deficient  (by  3),  and 
hence  there  will  be  infmitely-many  solutions  to  this  system. 
Successful  recovery  of  the  Fourier  phases  modulo  2 n  requires 
a  solution  preserving  the  integrality  of  the  error  term  .  In 
this  case,  we  obtain  a  final  RSC  solution  2 nc  away  from  the 
true  solution  for  some  integer  vector  c.  To  achieve  this,  we 
rely  on  the  integer  matrix  decomposition  known  as  the  Smith 
Normal  Form ,  which  is  described  in  the  following  Theorem 
and  Lemma: 

Theorem  2.2  (Smith  Normal  Form)  (Smith  1861):  Let  A 
be  a  nonzero  m  X  n  integer  matrix  with  rank  r.  There  ex¬ 
ist  integer  and  unimodular  4  (and  thus  invertible)  matrices 


3  Without  loss  of  generality,  we  have  selected  =  0  in  Equation 
(9)  so  as  to  simplify  the  diagram 

4  A  unimodular  matrix  is  one  whose  determinant  is  ±1 
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Noiseless 


Noisy 


/?  (wrapped) 


Figure  5.  Illustration  of  the  fundamental  ambiguity  of  27r-periodicity  in  RSC  imaging.  Distinct  unwrappings  /?*  and  p J  both  produce 
solutions  to  Equation  (4)  in  the  noiseless  case,  as  do  their  respective  projections  onto  K  +  L,  P*K+L  and  k+l  >  *n  the  n°isy  case. 


m  X  m  and  n  X  n  matrices  U  and  V  respectively  such  that 
the  matrix  product  D  =  UAV  is  a  diagonal  matrix  whose 
diagonal  entries  (the  so-called  elementary  divisors )  are 
non-zero  integers  for  i  <  r,  and  zero  for  i  >  r.  Moreover,  the 
matrices  U  {and  respectively,  the  matrix  V}  represent  the 
row  {and  column}  operations  which  zeroize  A  below  {and 
above}  the  diagonal.  □ 

Lemma  2.3  (Elementary  Divisors)  (Smith  1861): 
The  product  of  the  elementary  divisors  is  the  great¬ 
est  common  divisor  (gcd)  of  all  r  Xr  minors  of  A. 

The  proof  of  the  Theorem  and  Lemma  above  can  be 
found  in  textbooks  such  as  Newman  (1972).  □ 

Let  us  compute  the  Smith  Normal  Form  (SNF)  {U,  D,  V} 
of  our  matrix  M.  Let  U m  -  U_1  and  Mm  =  V-1  so  that  we 
can  write: 

M  =  U  mDmVm  (14) 

where  the  r  diagonal  elements  {dp  of  Dm  are  the  elementary 
divisors  of  M. 

We  can  now  re-write  Equation  (13)  above  as: 


aperture  pattern  are  all  1,  the  RSC  solution  defined  by  trsc 
is  immune  to  phase- wrapping  error. 

□ 

RSC  patterns  consisting  of  apertures  placed  randomly 
on  a  Cartesian  grid  appear  to  satisfy  this  sufficient  condition 
with  high  probability.  We  conducted  a  simple  experiment 
in  which  15  apertures  were  randomly  placed  on  a  10  X  10 
grid.  Out  of  256  placements,  66  were  valid  RSC  patterns 
(i.e.  possessed  at  least  critical  redundancy),  and  of  these, 
only  2  had  non- unity  elementary  divisors. 


2.3  Relation  to  closure-phase  approaches 

The  SNF  has  been  been  applied  to  the  RSC  phase  prob¬ 
lem  before  (Lannes  2003).  Whereas  we  have  chosen  to  ap¬ 
ply  SNF  directly  to  the  baseline  measurement  matrix,  the 
approach  taken  by  Lannes  (2003)  is  to  instead  treat  the 
piston- invariant  phases  of  the  bispectra  (the  so-called  clo¬ 
sure  phases)  as  the  fundamental  observables  from  which  the 
object  phases  can  be  inferred  via  the  relation: 


DmM  M 


e 

<f> 


-  u  mP*k+l 


(15) 


Let  us  choose  the  following  particular  solution  to  Equa¬ 
tion  (15): 

TRSC  =  VMDM^M^+L  (16) 

where  denotes  the  pseudo-inverse  of  D.  The  resulting 
error  is  then: 


eRSc  =  ^D+MU^e*h  (17) 

Lemma  2.4:  Let  u  =  e*  .  The  residual  wrap  error  epsc 

equals  0  mod  2 n  if  and  only  if  mod  {u^dp  =  0,  V/  <  r.  The 
proof  of  this  Lemma  is  an  adaptation  of  a  standard  proof 
which  can  be  found  in  most  textbooks  covering  linear  Dio- 
phantine  equations  (see,  e.g.  Newman  (1972)).  □ 

From  this  Lemma,  the  following  Corollary  is  clear: 
Corollary  2.5  (Sufficient  condition  on  SNF  of  RSC 
matrix  for  wrap -invariance):  If  the  elementary  divisors 
of  the  measurement  matrix  M  corresponding  to  a  certain 


Co^cO  =  yci  +2tt  ect  (18) 

where  yci  are  the  wrapped  closure  phases,  ec;  is  the  corre¬ 
sponding  wrap  vector,  and  C0->c  is  the  matrix  mapping  the 
distinct  object  phases  in  the  array  to  closure  phases.  Lannes 
(2003)  hence  applies  the  SNF  to  the  closure  matrix  C0^c.  By 
direct  analogy  to  Corollary  2.5,  note  that  if  the  elementary 
divisors  of  Coc  are  all  1,  then  the  RSC  solution  is  immune 
to  phase- wrapping  error.  Using  closure  phases  as  observables 
can  be  advantageous  in  low-light  scenarios  in  which  there  is 
not  sufficient  SNR  in  a  single  atmospheric  coherence  time 
to  reliably  measure  the  baseline  phases.  To  overcome  this 
low  per-frame  SNR,  atmosphere- invariant  observables  such 
as  the  bispectra  can  be  integrated  over  many  frames  to  build 
sufficient  SNR,  and  their  respective  closure  phases  used  as 
reliable  phase  measurements.  Since  the  baseline  phases  are 
known  modulo  2n,  the  linear  combinations  of  them  which 
comprise  the  closure  phases  are  also  known  modulo  2/r.  In 
order  to  relate  this  condition  to  Corollary  2.5,  let  us  first 
define  another  closure  matrix  Cm_>c  which  instead  maps  the 
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Figure  6.  Distinction  between  spanning  tree  baselines  (thick, 
solid)  and  loop  entry  baselines  (thin,  dotted) 


phase  measurements  to  closure  phases.  This  mapping  con¬ 
sists  of  equations  of  the  form: 


y  123  =  Pl2  +  023  “  P\3  (19) 

where  y 123  is  the  closure  phase  associated  with  apertures 
1,  2,  and  3,  and  the  fiij  are  the  associated  baseline  phases 
(see  Equation  (1)).  Of  the  possible  closure  phases,  at 

most  can  be  linearly-independent  (see  e.g.  Readhead 

et  al.  (1988)).  One  commonly-chosen  set  of  such  linearly- 
independent  relations  consists  of  all  the  closure  triangles  in¬ 
volving  a  given  aperture  A,  and  this  is  the  set  selected  by 
Lannes  (2003).  Cm_>c  is  therefore  an  x  (2O  matrix- 

Lannes  (2003)  accordingly  provides  a  convenient  grouping 
of  the  baselines  into  two  categories:  (1)  spanning  tree  base¬ 
lines  which  connect  aperture  A  to  all  other  apertures,  and 
(2)  loop  entry  baselines  which  provide  the  closure  for  these 
spanning  tree  baselines.  This  categorization  is  depicted  in 
Figure  6. 

Given  this  categorization,  we  can  decompose  Cm_»c  into 
corresponding  blocks  as: 

c m->c  =  [  Cm_>c  I  (AT-1)  |  (20) 

where  Cm_>c  contains  the  spanning  tree  contributions  to  the 
matrix  (which  appear  in  multiple  closures),  and  \,n-i\  is  the 

X  identity  matrix  representing  the  loop-entry 

contributions  (each  of  which  appears  in  only  one  closure). 
The  following  property  follows  from  this  block  form  expres¬ 
sion: 

Lemma  2.6:  The  elementary  divisors  of  Cm_>c  are  all  1. 
Proof:  Since  we  have  chosen  a  linearly-independent  subset 
of  closure  relations,  r  =  rank(Cm^c)  =  (^2^)-  There  exists 
a  r  Xr  minor  (namely,  Ljv-k)  which  is  equal  to  1.  Therefore 
the  gcd  of  all  r  X  r  minors  is  1,  and  therefore  from  Lemma 
2.3  (Elementary  Divisors),  all  elementary  divisors  must  be 
1.  □ 
Let  us  now  relate  Cm^c  to  the  matrix  C0^c  used  by 
Lannes  (2003).  Recall  from  the  discussion  of  bispectra  in  Sec¬ 
tion  1  that  the  closure  relations  eliminate  piston  differences 
in  the  measurements  so  that  Cm_>c  annihilates  the  subspace 
L,  i.e.  the  space  spanned  by  the  columns  of  M  correspond¬ 
ing  to  0.  Defining  M#  as  the  submatrix  of  M  containing  the 
columns  corresponding  to  6 ,  we  have 


where  C0_>c  =  C0_>c  is  an  ^xd  matrix  which 

is  rank-deficient  by  two  5. 

Then  by  direct  analogy  to  Equation  (13),  we  can  obtain 
valid  RSC  object  phase  solutions  by  solving: 


e 

0 


[  >c  0  j 


e 

0 


-  yci +  2 neh,ci  (22) 


where  y*cl  is  the  true  unwrapped  closure  vector  and 
2;re/z  cl  *s  residual  integer  wrapping  error  vector  after 
applying  the  Babai  NP  algorithm  to  solve  the  CVP  prob¬ 
lem  associated  with  matrix  Cm^cMg  and  2/rec/.  Note  that 
if  we  find  a  vector  6 *  satisfying  Equation  (22),  it  will  clearly 
also  satisfy  the  relation  Co_>c0*  =  y*cl  +  2 ne^  cl  of  Lannes 
(2003).  Note  furthermore  that  we  can  solve  the  equation 
above  in  two  separate  integer-preserving  steps  of  the  form 
of  Equation  (16),  the  first  involving  the  SNF  decomposition 
of  Cm_>c,  and  the  second  involving  that  of  M.  Since  the 
elementary  divisors  of  Cm_>c  are  all  1  by  construction  (by 
Lemma  2.6),  and  hence  the  first  step  is  integer-preserving, 
wrap- invariance  can  again  be  ensured  if  all  the  elementary 
divisors  of  M  are  1.  Therefore  we  have  the  following  Propo¬ 
sition  for  wrap- invariance  in  closure-based  estimation  of  the 
Fourier  phase: 

Proposition  2.7  (Sufficient  condition  for  wrap- 
invariance  of  closure-based  RSC):  If  the  elementary  di¬ 
visors  of  the  phase  measurement  matrix  M  are  all  1,  then  the 
closure-based  RSC  solution  will  be  wrap- invariant.  Proof: 
(see  Appendix  A)  □ 

We  remark  in  passing  that  although  the  preceding  anal¬ 
ysis  was  presented  in  the  context  of  the  traditional  three- 
aperture  closure,  it  applies  directly  to  the  case  of  closures 
involving  an  arbitrary  number  of  sides.  As  an  example,  con¬ 
sider  the  pattern  shown  in  Figure  7.  A  spanning  tree  for 
the  pattern  consisting  of  the  short  baselines  in  an  array  is 
depicted.  Let  {(f)sp}  denote  the  aperture  phase  differences  in 
these  N  -  1  spanning  tree  baselines.  Note  that  all  aperture 
phase  differences  in  the  array  can  be  expressed  as  linear  com¬ 
binations  of  the  {(psP)-  H  the  aperture  phase  differences  are 
known  reliably  via  measurements  of  the  spanning  tree  base¬ 
lines,  we  can  use  these  measurements  to  cancel  the  aperture 
phase  differences  in  all  other  measurements  (of  which  one  ex¬ 
ample  is  shown  in  green).  The  idea  of  using  such  generalized 
closure  phases  (Martinache  2010)  is  indeed  the  mathemati¬ 
cal  foundation  for  the  promising  technique  known  as  baseline 
bootstrapping  in  which  high-fidelity  phase  measurements  of 
several  high- SNR  baselines  (typically  the  short  baselines)  to 
cancel  the  atmosphere  on  each  lower- SNR  (and  hence  lower 
fidelity)  baseline. 

Note  that  for  an  arbitrary  A-aperture  pattern,  there 
will  in  general  be  N  -  1  spanning  tree  baselines  and  thus 
( 2O  “  (N  -  1)  =  (A2_1)  generalized  closures,  each  involving 
a  distinct  closing  (or  loop-entry)  baseline.  Therefore  the  re¬ 
sulting  measurement  matrix  can  be  expressed  exactly  as  in 
Equation  (20)  above  and  hence  the  preceding  analysis  holds. 

While  this  section  has  considered  a  few  possibilities 
for  phase  observables,  relating  mathematical  conditions  for 


6 

0 


(21) 


5  To  see  this,  note  that  each  solution  set  to  Equation  (18)  above 
remains  valid  after  replacing  each  Oij  with  0 f.  =  dij  -  z  •  (rj  -  rp 
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Figure  7.  Bootstrapping  phase  of  a  low-SNR  baseline  (green) 
with  subset  (blue)  of  high-SNR  baselines  from  spanning  tree  base¬ 
lines  (black) 

wrap  invariance  to  a  physical  condition  on  aperture  place¬ 
ment  is  more  intuitive  when  considering  the  raw  phase  mea¬ 
surements  as  opposed  to  their  closures.  For  this  reason,  for 
the  remainder  of  the  paper  we  will  work  directly  with  the 
baseline  phases  and  their  associated  wrapping  errors.  In  par¬ 
ticular,  we  will  begin  by  connecting  these  wrapping  errors 
with  analogous  ambiguities  in  recently-developed  phasor- 
based  approaches. 


3  IMPLICATIONS  OF  WRAP  AMBIGUITIES 

ON  PATTERN  DESIGN 

In  this  section  we  use  the  mathematically-sufficient  condi¬ 
tions  for  wrap  invariance  from  the  previous  section  to  show 
that  aperture  patterns  whose  interferometric  graph  satisfies 
a  certain  loop- free  condition  are  wrap- invariant.  Here  we 
define  the  interferometric  graph  in  the  standard  way:  it  is 
simply  the  graph  formed  by  connecting  the  array’s  apertures 
(the  nodes  of  the  graph)  with  edges  representing  the  array’s 
baselines. 

This  condition  is  founded  on  the  Lemma  2.3  (Elemen¬ 
tary  Divisors)  in  Section  2  and  the  following  definition  of  the 
matrix  determinant.  This  definition  is  given  in  many  linear 
algebra  texts  (see  e.g.  Bretscher  (2001)). 

Definition  3.1  :  Suppose  we  have  an  nxn  matrix  A.  Define 
a  pattern  as  a  selection  of  n  entries  of  the  matrix  in  which 
there  is  only  one  chosen  entry  in  each  row  and  one  in  each 
column  of  the  matrix.  Furthermore,  we  denote  a  pair  of  num¬ 
bers  in  a  pattern  as  inverted  if  one  of  them  is  located  above 
and  to  the  right  of  the  other.  Then  we  can  obtain  the  deter¬ 
minant  of  A  by  summing  the  products  associated  with  all 
patterns  with  an  even  number  of  inversions  and  subtracting 
the  products  associated  with  all  the  patterns  with  an  odd 
number  of  inversions. 

Consider  a  r  X  r  sub-matrix  M /  consisting  of  a  set  I 
of  linearly-independent  rows  in  M  and  d  +  N  -  3  linearly- 
independent  columns.  Our  goal  will  be  to  find  conditions 
under  which  such  a  sub-matrix  contains  only  one  pattern 
with  a  non-zero  product,  in  which  case  the  determinant  will 
be  +1  from  the  definition  above.  This  will  guarantee,  by  the 
Lemma  2.3  (Elementary  Divisors)  of  the  previous  section, 
that  the  elementary  divisors  of  M  are  all  1,  which  will  in 
turn  ensure  that  the  error  in  the  RSC  solution  trsc  will  be 
0  mod  2 n  by  Corollary  2.5. 

We  begin  our  search  for  such  a  matrix  by  forming  the 


intermediate  r  X  (r  +  3)  matrix  M /  consisting  of  I  linearly- 
independent  rows  in  M,  and  define  the  corresponding  sub¬ 
graph  G  containing  only  the  baselines  corresponding  to  these 
rows.  In  particular,  let  us  examine  the  pattern  restrictions 
encountered  if  we  select  a  given  column  in  M  /  for  partici¬ 
pation  in  our  r  X  r  sub- matrix  M  /.  Namely,  we  will  sequen¬ 
tially  identify  those  special  matrix  entries  which  must  be 
part  of  a  pattern  with  a  non-zero  product.  For  example,  all 
non-redundant  measurements  contain  a  singleton  +1  in  the 
column  associated  with  their  object  phase.  All  non-zero  pat¬ 
terns  must  clearly  contain  this  +1  and  so  we  can  select  these 
singleton  object-phase  entries  as  guaranteed  participants  in 
a  non-zero  pattern.  Moreover,  all  measurements  containing 
a  leaf  node  (i.e.  a  node  with  a  single  connection)  in  G  contain 
a  singleton  +1  in  the  column  associated  with  their  leaf  node. 
All  non-zero  patterns  must  clearly  contain  this  +1  as  well. 
Thus  we  can  also  select  these  leaf  node  entries  as  guaranteed 
participants  in  a  non-zero  pattern. 

There  may  be  cascading  implications  of  such  single- 
ton  measurements.  To  illustrate  this,  consider  the  scenario 
shown  in  Figure  8.  A  simple  RSC  array  is  shown  on  the 
left.  A  subset  of  the  baselines  in  one  possible  linearly- 
independent  sub-matrix  M/  is  depicted.  A  simplified  depic¬ 
tion  of  the  matrix  M/  is  shown  in  which  all  non-zero  entries 
have  been  colored  black  and  all  zero  entries  have  been  col¬ 
ored  white  for  simplicity.  In  Step  A  of  the  reduction  process, 
object  phase  65  is  selected  for  participation  since  it  is  a  sin¬ 
gleton  object  phase.  Its  corresponding  row  (i.e.  row  5)  in 
M  /  is  then  eliminated  from  participation  since  the  remain¬ 
ing  entries  in  this  row  cannot  participate  in  a  pattern  (by 
definition  of  a  pattern).  In  Step  B,  the  entry  correspond¬ 
ing  to  aperture  6  (i.e.  06 )  in  row  9  is  selected  since  it  has 
become  a  leaf  node  in  the  pattern,  and  row  9  can  then  be 
eliminated.  Then  in  Step  C,  object  phase  64  is  then  selected 
by  virtue  of  becoming  a  singleton  object  phase,  and  row  4  is 
then  eliminated.  This  selection/elimination  process  can  be 
repeated  beyond  the  steps  shown  in  the  Figure,  until  either 
no  leaf  nodes  and  singleton  object  phases  remain,  or  there 
are  no  baselines  left  to  eliminate.  We  formalize  the  pattern 
reduction  process  in  Algorithm  1  below. 


Algorithm  1  Pattern  Reduction  Algorithm 
Require:  R  {where  R  is  the  set  of  the  baselines  corre¬ 
sponding  to  M/,  where  I  denotes  the  indices  of  a  linearly- 
independent  subset  of  d  +  N  -  3  rows  of  M } 

while  \R\  >  0  do 

1.  Leaf  Nodes 

1.1  remove  any  remaining  baselines  containing  leaf 
nodes  from  R 

1.2  add  the  associated  apertures  to  the  list  N 

2.  Singleton  Object  Phases 

2.1  remove  any  remaining  baselines  containing  singleton 
object  phases  from  R 

2.2  add  the  associated  object  phases  to  the  list  O 

if  no  baselines  removed  in  the  current  iteration  then 
return  PERSISTENT 
end  if 
end  while 
return  LOOPFREE 


We  can  see  that  any  baseline  in  an  interferometric  graph 
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Figure  8.  Example:  Reducing  an  aperture  pattern  and  associated  matrix  to  identify  Persistent  Loop(s) 


which  does  not  belong  to  a  loop  will  be  eliminated  in  the 
reduction  process,  and  its  corresponding  matrix  entries  fac¬ 
tored.  The  only  structures  in  the  graph  that  persist  after  this 
reduction  are  sets  of  loops  with  a  certain  property.  Namely 
we  define  a  persistent  loop  set  as  a  set  of  loops  that  contains 
at  least  two  instances  of  every  baseline  contained  in  the  set. 
(A  set  can  consist  of  any  number  of  loops,  including  one). 
With  this  definition,  absolute  invariance  may  be  possible  if 
the  graph  of  the  redundant  baselines  does  not  contain  any 
persistent  loop  sets.  Algorithm  1  returns  PERSISTENT  if 
persistent  loops  exist  and  LOOPFREE  if  the  pattern  is  com¬ 
pletely  reduced  and  therefore  free  of  persistent  loops. 

Note  that  in  the  latter  case,  we  will  have  eliminated 
r  rows  from  M7.  Since  each  baseline  elimination  is  associ¬ 
ated  with  object  phase  or  aperture  selection  from  distinct 
columns,  and  M/  contains  r  +  3  columns,  there  will  be  ex¬ 
actly  three  extraneous  columns  not  involved  in  the  reduction 
process.  The  r  Xr  submatrix  obtained  by  selecting  the  non- 
extraneous  columns  (i.e.  those  corresponding  to  the  selected 
object  phases  and  leaf  nodes  in  O  and  A,  respectively)  will 
then  have  a  determinant  of  +1  by  virtue  of  having  a  single 
non-zero  pattern  revealed  by  the  reduction  process.  Having 
ensured  the  existence  of  a  unit  r  X  r  minor,  the  gcd  of  all 
r  x  r  minors  must  be  1.  We  have  hence  confirmed  the  ele¬ 
mentary  divisors  must  be  all  1,  and  thereby  ensured  that 
the  pattern  is  wrap- invariant  (c.f.  Lemma  2.3  (Elementary 
Divisors),  and  Corollary  2.5).  We  summarize  the  sufficient 
condition  as  follows: 

Proposition  3.2  (Sufficient  conditions  on  aperture 
pattern  for  wrap -invariance):  Consider  the  graph  of  an 
aperture  pattern  which  contains  d  distinct  baselines  and  any 
set  of  N  -  3  linearly- independent  redundant  baselines.  If  this 
graph  does  not  contain  persistent  loop  sets  (in  the  sense  de¬ 
fined  above),  the  matrix  M/  formed  by  these  independent 
measurements  will  have  determinant  +1.  As  a  result  Corol¬ 
lary  2.5  will  hold,  thereby  guaranteeing  that  RSC  solution 


Persistent  Loop 


Figure  9.  Persistent  Loop  set  at  center  of  pattern  in  Figure  2 

trsc  wiH  be  invariant  to  wrapping  of  the  phase  measure¬ 
ments. 

We  have  hence  arrived  at  a  a  physical  definition  of  a 
wrap-invariant  pattern.  We  now  apply  Algorithm  1  to  the 
example  pattern  shown  in  Figure  2.  Algorithm  1  reduces  the 
pattern  to  the  persistent  loop  set  shown  in  Figure  9. 

The  elementary  divisors  of  the  pattern’s  measurement 
matrix  are  not  all  1;  they  are  all  1  except  for  a  singleton  3 
and  hence  mod  3  =  0  for  all  choices  of  M/. 

Having  traced  the  distortion  induced  by  phase  wrapping 
to  physical  property  of  the  array  itself,  we  now  return  to 
Figure  3.  The  lower  left  panel  shows  the  reconstruction  result 
with  the  SNF-based  Phase  method  described  in  Section  2. 
The  closure  phase  approach  yields  the  same  corruption  in 
reconstruction,  as  the  elementary  divisors  of  Coc  are  also  all 
1  except  for  a  singleton  3. 

There  are  several  simple  ways  to  amend  this  pattern 
so  that  it  is  wrap- invariant.  While  the  most  intuitive  of 
these  involve  moving  the  apertures  involved  in  the  persis¬ 
tent  loop  shown  in  Figure  9,  these  approaches  leave  gaps  in 
the  UV-sampling  pattern.  An  alternate  approach  that  pre¬ 
serves  the  UV-sampling  is  to  add  an  aperture  to  the  center  of 


MNRAS  000,  1-15  (2015) 


10  B.  G.  Kurien  et  al. 


Figure  10.  Amended  Pattern 
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Figure  11.  Reconstruction  Results  for  Amended  Pattern 
Phase  Approach  (Noiseless)  Phase  Approach  (SNR  =  25  dB) 


the  pattern  as  shown  in  Figure  10.  This  results  in  additional 
linearly-independent  redundancies  colored  in  blue  and  green, 
respectively,  in  the  Figure.  These  additions  replace  baselines 
in  the  persistent  loop,  allowing  this  loop  to  be  broken.  With 
wrap-invariance,  reconstruction  results  match  the  true  im¬ 
age  in  both  the  phase  and  phasor  approaches  as  respectively 
shown  in  Figure  11.  In  the  top  row,  reconstruction  results  are 
displayed  for  the  phase  (left)  and  phasor  (right)  approaches 
for  the  noiseless  case.  The  image  distortion  present  in  Figure 
3  has  been  completely  eliminated  by  tweaking  the  pattern 
so  that  it  is  wrap- invariant.  Analogous  results  for  an  SNR  of 
25  dB  are  displayed  in  the  bottom  row.  Here  we  define  SNR 
as  the  ratio  of  the  phasor  magnitude  at  visibility  1  (i.e.  zero 
spatial  frequency)  to  the  standard  deviation  of  the  noise, 
which  we  have  assumed  to  be  complex  Gaussian  and  i.i.d. 
across  spatial  frequency  for  this  simulation. 


4  WRAP-INVARIANCE  AND  PRACTICAL 
RSC  CALIBRATION 

In  previous  sections  we  established  that  the  Phase  Approach 
can  be  made  robust  to  phase- wrapping  using  the  Smith  Nor¬ 
mal  Form  (SNF)  and  algorithms  from  lattice  theory.  More¬ 
over,  the  SNF  provided  a  mathematical  framework  for  the 
notion  of  wrap- invariance.  From  a  practical  standpoint,  how¬ 
ever,  computation  of  the  Smith  Normal  Form  is  likely  to  be¬ 
come  a  computational  burden  for  large  arrays,  as  in  those 
under  current  consideration  with  Nap  «  102  and  n  «  104 
baselines  (Zheng  et  al.  2014).  Techniques  not  requiring  such 


a  computation  are  hence  of  strong  practical  interest.  In  this 
section  we  show  that  the  wrap- invariance  property  checked 
by  Algorithm  1  provides  a  certificate  for  reliable  reconstruc¬ 
tion  with  such  techniques,  in  the  presence  of  the  fundamental 
2/r-periodicity  of  interferometric  measurements. 

4.1  Practical  Phase  Approaches 

Our  approach  here  will  rely  upon  the  well-known  Singular 
Value  Decomposition  (SVD)  of  the  measurement  matrix  M, 
which  is  given  by: 

M  =  VcrLrVl  (23) 

in  which  Ucr  and  V^-  are  m  X  m  and  ( d  +  N)  X  (d  +  N)  or¬ 
thogonal  matrices,  respectively.  Lo-  is  a  m  x  ( d  +  N)  diagonal 
matrix  with  r  non-zero  diagonal  entries  (the  so-called  singu¬ 
lar  values  of  M),  where  r  =  rank(M)  =  d  +  N  -  3. 

Lemma  4.1:  The  final  3  columns  of  V^-  form  a  basis  for  the 
nullspace  of  M. 

Proof:  This  follows  from  the  fact  M  is  rank-deficient  by  3, 
and  standard  properties  of  the  right  singular  vectors  com¬ 
prising  V q-  in  the  SVD.  (Bretscher  2001)  □ 

Now  recall  that  in  Section  2,  we  provided  an  one  par¬ 
ticular,  SNF-based,  solution  to  Equation  (13).  Here  we  in¬ 
stead  consider  the  complete  set  of  solutions  to  this  Equation, 
which  is  given  by: 

Jo-  —  M  (/^o ,K+L  +  T0  (24) 

where  tq  is  any  vector  in  the  nullspace  of  M,  and  M+  denotes 
the  pseudo-inverse  of  M,  whose  matrix  elements  are  derived 
directly  from  the  SVD  above. 

M+=V(TE+U^  (25) 

where  is  a  ( d  +  N)  X  m  diagonal  matrix  whose  r  non-zero 
diagonal  entries  are  the  reciprocals  of  the  corresponding  non¬ 
zero  entries  in  Lo-. 

Typically  Phase  Approach  techniques  implicitly  select 
a  particular  solution  from  the  family  of  solutions  in  Equa- 
ton  (24)  via  augmenting  the  matrix  M  with  additional  con¬ 
straints.  The  most  common  among  these  (Wieringa  1992) 
(Wijnholds  &  Noorishad  2012)  is  to  enforce  that:  2  0  =  0, 
2  Vi  =  01  where  r i  is  the  vector  position  of  the  i-th  aper¬ 
ture  in  the  array. 

Note  that  the  error  resulting  from  application  of  this 
pseudo-inverse  to  the  unwrapped  measurement  vector  will 
be  given  by: 

2/recr  =  M+(2^ep  (26) 

Note  furthermore  that  the  vector  representing  the  solu¬ 
tion  error  ea-  has  two  parts:  one  for  the  error  in  the  Fourier 
phases  (which  we  denote  with  the  index  K),  and  one  for  the 
error  in  the  atmosphere/piston  (which  we  will  denote  with 
the  index  L),  i.e.: 

£ct  =  [ eo-,K ,  Co-,l]T  (27) 

We  will  focus  attention  on  the  latter,  since  it  is  of  direct 
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relevance  for  image  formation.  Let  us  express  the  spatial 
frequencies  measured  by  an  array  as  two-element  vectors  of 
the  form  (ajx,cjjy).  Let  X  be  the  dx 2  matrix  containing  these 
spatial  frequencies.  Note  then  that  the  phase- wrap  error  will 
manifest  itself  merely  as  an  image  shift  if  and  only  if  this 
error  is  a  (modulo-2^-)  phase  ramp,  i.e.  there  exists  a  2- 
element  shift  vector  z  and  an  integer  vector  k  which  together 
satisfy: 

2n e<j-,K  ~  2 nXz  =  Ink  (28) 

Substituting  from  Equation  (26)  we  obtain: 

M^(2/re)  -  2 nXz  =  Ink  (29) 

where  denotes  the  sub-matrix  of  M+  formed  by  the  rows 
associated  with  K. 

Dividing  through  by  2 n  we  obtain  the  equation:  M^e* - 
Xz  =  k.  Note  that  each  element  of  can  be  expressed 
as  some  rational  number  Similarly  we  first  assume  X 
contains  rational  spatial  frequencies  with  greatest  common 
denominator  qx .  Then  we  can  multiply  through  by  the  least- 
common- multiple  (LCM)  of  the  { qi }  and  qx  to  obtain  a  sys¬ 
tem  of  equations  whose  coefficients  are  guaranteed  to  be 
integer  (i.e.,  we  have  a  linear  Diophantine  system).  Let  this 
LCM  be  denoted  as  /.  Then  we  have,  after  rearranging  terms, 

IXz  =  /(M+e*  -fc)  (30) 

We  now  wish  to  determine  conditions  under  which  there 
exist  vectors  k  and  z  satisfying  this  overdetermined  Dio¬ 
phantine  system.  Applying  the  Smith  Normal  Form  decom¬ 
position  (c.f.  Theorem  2.2)  to  the  matrix  IX  this  time,  and 
noting  that  rankiX )  -  2,  we  have: 

DX  =  Ux(/X)Vx  (31) 

where  ux  and  vx  are  unimodular  matrices  of  size  dxd 
and  2x2,  respectively,  and  Dx  is  a  rectangular  diagonal 
matrix  whose  entries  are  zero  below  row  2. 

If  we  left-multiply  Equation  (30)  by  Ux  on  both  sides, 
we  obtain: 

/UxX^  =  /Ux(M+e*  -k)  (32) 

Using  Equation  (31)  and  the  fact  that  Vx  is  a  unimod¬ 
ular  (and  hence  invertible)  matrix,  we  can  then  write: 

DxVx1z  =  /(UxM+e;-Uxfe)  (33) 

We  are  now  in  position  to  prove  the  main  result  of  this 
section,  which  is  preceded  by  the  following  Lemma: 
Lemma  4.2:  Given  wrap- invariance,  the  (column)  vector 
UXM+ e*  has  integer  entries  below  row  2. 

Proof :  (see  Appendix  B)  □ 

Proposition  4.3:  If  a  pattern  is  wrap- invariant  (in  the  sense 
of  Section  3),  reconstruction  error  induced  by  phase  wrap¬ 
ping  is  limited  to  an  image  shift. 

Proof :  We  re-arrange  the  Equation  (33)  above  so  that  it 
reads: 

yDxV^z-UxM+e^-Uxfc  (34) 


Let  v  =  jDxV^z  -  UxM^e^.  Note  that  since  DX  is 
zero  below  row  2,  the  entries  of  v  below  row  2  will  be  equal 
to  those  of  (-  UXM+e;  ),  which  are  integers  by  Lemma  4.2. 
Now  consider  the  first  and  second  entries  of  v.  Let  f  be 
the  vector  containing  the  fractional  parts  of  the  first  two 
elements  of  vector  UxM^e^,  and  let  A  be  the  invertible 
matrix  consisting  of  the  first  two  rows  of  j  Dx V^1 .  Without 
loss  of  generality,  choose  z*  =  A-1/  so  that  the  fractional 
part  f  is  annihilated,  leaving  only  integer  elements  in  the 
first  two  entries  of  v.  Hence  we  now  have: 

v  =  -Uxk  (35) 

with  v  ensured  to  contain  only  integer  elements.  Since 
Ux  is  unimodular,  the  vector  k*  =  -U^1^  will  be  integral. 
We  have  thus  found  a  pair  ( z*,k *)  with  integer  fc*  which 
satisfies  the  Equation  (33).  Since  Equation  (33)  is  related 
to  Equation  (30)  via  a  unimodular  (and  hence  invertible) 
mapping  Ux,  invariance  is  hence  proven.  □ 

With  the  previous  result,  we  have  characterized  the 
complete  family  of  Phase  Approach  solutions  given  in  Equa¬ 
tion  (24).  Namely  we  have  shown  that,  for  a  wrap- invariant 
pattern,  the  family  differs  by  at  most  an  image  shift  from  the 
true  solution.  Returning  to  our  running  example  in  Figure 
10,  we  verified  that  different  solution  choices  from  the  family 
given  in  Equation  (24)  simply  resulted  in  shifts  of  an  other¬ 
wise  pristine  image  in  the  reconstruction.  On  the  other  hand 
with  wrap- variant  pattern  in  Figure  2,  image  distortion  of 
the  severity  of  Figure  3  was  again  observed,  as  expected. 

4.2  Practical  Phasor  Approaches 

Though  traditional  treatments  employ  the  phase  approach 
of  the  previous  section  which  operates  on  baseline  phases, 
recent  papers  (e.g.  Marthi  &  Chengalur  (2014),  Liu  et  al. 
(2010))  have  shown  that  approaches  which  operate  at  the 
phasor  level  can  be  superior  in  accuracy.  Liu  et  al.  (2010) 
developed  a  Gauss-Newt  on- type  Non-linear  Least-Squares 
(NLS)  solver  and  showed  it  produced  unbiased  phase  es¬ 
timates,  in  contrast  with  the  biased  ones  provided  by  the 
phase  approach.  Marthi  &  Chengalur  (2014)  and  Wijn- 
holds  &  Noorishad  (2012)  have  also  proposed  low-complexity 
phasor-based  approaches  and  demonstrated  performance 
near  the  Cramer- Rao  Bound.  Though  the  capacity  of  the 
Phasor  approaches  to  produce  superior  accuracy  relative  to 
Phase  approach  has  been  demonstrated,  the  former’s  con¬ 
vergence  issues  can  be  mitigated  via  initialization  with  the 
results  of  the  latter  (see  e.g.  Liu  et  al.  (2010),  and  Zheng 
et  al.  (2014)). 

The  implementations  of  the  Phasor  Approach  typically 
employ  the  following  measurement  model: 

Vij  =  SiS*jfij  +  nij  (36) 

where  Uj  is  the  complex  visibility  observed  between  aper¬ 
tures  i  and  j,  gi  =  \gi\eJ^i  and  gj  =  \gj\eJ(f>J  are  the  complex 
gains  of  these  apertures,  fij  is  the  true  complex  visibility 
measured  by  this  pair,  and  n  is  complex  measurement  noise. 
Note  that  the  phase  difference  between  gi  and  gj  is  sim¬ 
ply  the  optical  path  difference  between  apertures  i  and  j 
introduced  in  the  previous  section.  Given  this  model,  NLS 
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approaches  attempt  to  find  a  set  of  complex  phasors  {g*}  and 
{/}  which  minimize  an  objective  function  of  the  form: 

A  =  TjTj  wiMVij  -  8i8jfij)(V*j  ~  g-gjfij)  II  (37) 

i  j>i 

Minimization  of  A  with  respect  to  the  unknowns  (i.e. 
distinct  object  and  antenna  complex  gains)  can  be  accom¬ 
plished  with  iterative  application  of  the  following  updates, 
as  reported  by  Marthi  V  Chengalur  (2014)  in  the  context 
of  radio  interferometry,  and  by  Lacour  et  al.  (2007)  in  the 
context  of  optical  interferometry: 


2 j ±k  wkj  gj  fkj-  Vkj 

8k  ~  ljj*kwkj\gj\2\fkj\2 

(38) 

_  Yuj>kgfrgjVkj 
b~  Zj>kWkj\gk\2\gj\2 

(39) 

where  the  {/^}  are  the  true  complex  visibilities  of  the  distinct 
object  phases  in  the  array. 

Due  to  the  circularity  of  these  definitions,  these  equa¬ 
tions  must  be  solved  iteratively.  Starting  from  an  initial 
guess  for  all  phasors,  Equation  (38)  is  solved  to  obtain  a 
better  estimate  for  the  {g^}  and  then  these  {g^}  are  used  to 
obtain  refined  estimates  of  the  {/^}  through  Equation  (39). 
In  the  next  iteration,  these  {fjy}  are  used  to  further  refine 
{g^},  and  so  on. 

Though  there  are  other  means  of  minimizing  objectives 
of  the  form  A  (Wijnholds  &  Noorishad  2012)  (Liu  et  al. 
2010),  we  omit  discussion  of  them  here;  our  present  purpose 
is  to  characterize  the  correctness  of  the  solutions  themselves, 
regardless  of  how  they  are  obtained.  As  has  been  noted  be¬ 
fore  (see  e.g.  Lannes  V  Anterrieu  (1999)),  there  are  strong 
connections  between  phase-  and  phasor-based  approaches. 
To  see  this,  let  z  be  the  vector  of  products  {gig*-  f\i-j\}  which 
minimize  A.  We  rewrite  Equation  (37)  as: 

A  =  E  2)  "y  II Wj  -  zij)(v*j  -  z*j) II  (40) 

i  j>i 

Define  rij  -  ej2nniJ  Zij  for  an  arbitrary  integer  nij  and  r 
as  the  vector  containing  the  r^j. 

Note  that  r  also  minimizes  A  since  the  rotations 
{ejlnriij  j  noj-  change  the  values  of  the  residuals  in  A. 
Hence  any  set  of  rotated  phasors  {g*}  and  {/|j_y|}  whose 
products  produce  the  vector  r  will  also  minimize  A.  Note 
that  the  set  of  such  valid  phase  vectors  (i.e.  the  concatena¬ 
tions  of  possible  {Zg/}  and  {Zf\i_j\})  includes  the  complete 
family  of  Phase- approach  solutions  To-  in  Section  4.1  with 
Pk+l  =  Zr  (where  Zr  is  the  vector  of  the  phases  of  the 
complex  vector  r).  In  other  words,  the  valid  phase  compo¬ 
nent  of  the  phasor  approach  solutions  is  not  unique,  and 
the  minimization  of  A  admits  the  same  solution  ambiguity 
depicted  in  Figure  5.  Hence  we  see  that  integer  ambiguities 
present  in  the  phase  approach  do  not  disappear  in  the  pha¬ 
sor  approach;  in  fact,  the  unwrapped  candidate  solutions 
of  the  phase-based  approach  correspond  to  minima  of  the 
phasor-based  objective.  In  practice  Phasor  Approach  tech¬ 
niques  may  converge  to  any  one  of  these  minima.  Hence  the 
critical  issue  for  reliable  image  reconstruction  is  again  the 


Figure  12.  Reconstruction  Results  for  Phasor  Approach 
Phasor  Approach  (Noiseless)  Phasor  Approach  (SNR  =  25  dB) 
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nature  of  the  difference  between  these  valid  minima  and  the 
true  solutions.  Based  on  the  connections  we  have  drawn  with 
Phase-approach  solutions  above,  the  following  Proposition  is 
clear: 

Proposition  4.4:  If  a  pattern  is  wrap- invariant,  the  global 
minima  of  the  Phasor- Approach  objective  caused  by  the  in¬ 
herent  2^-periodicity  in  the  objective’s  residual  differ  from 
the  true  solution  merely  by  an  image  shift. 

In  practice  we  see  that,  as  in  the  Phase  approach,  if  the 
pattern  is  not  wrap-invariant,  the  Phasor  Approach  suffers 
from  false  global  minima  producing  severe  distortion  in  the 
resulting  reconstruction.  The  lower  right  panel  of  Figure  3 
shows  the  resulting  reconstruction  produced  by  the  Phasor 
method  using  the  updates  in  Equations  (38)  and  (39).  To 
show  the  correspondence  in  solutions  between  two  types  of 
approaches,  we  provided  the  result  of  the  Phase  approach 
as  the  initial  point  for  the  Phasor  approach  as  is  common  in 
practice  (Liu  et  al.  2010),  (Zheng  et  al.  2014).  Indeed  this 
point  is  a  global  minimum  of  A  which  the  updates  above 
cannot  escape,  and  as  a  result  we  observe  virtually-identical 
distortion  to  that  of  the  Phase  approach. 

We  repeated  the  experiment  with  the  wrap- invariant 
pattern  in  Figure  10,  and  as  expected,  the  results  were  pris¬ 
tine  in  the  noiseless  case  and  virtually-identical  to  those  of 
the  Phase  approach  in  the  noisy  case.  For  completeness  the 
results  are  given  in  Figure  12. 

As  expected,  different  initializations  of  the  updates  in 
Equations  (38)  and  (39)  simply  resulted  in  shifts  of  an  oth¬ 
erwise  pristine  image  in  the  reconstruction. 


5  CONCLUSIONS 

In  this  paper,  we  have  examined  the  ambiguities  caused 
by  the  2/r-periodicity  of  interferometric  phase  in  Redundant 
Spacing  Calibration.  In  particular  we  have  described  their 
fundamental  presence  in  existing  RSC  methods  whether  the 
observables  considered  are  the  measured  baseline  phasors 
or  their  phases.  In  the  former,  e.g.  by  Greenaway  (1990), 
they  are  manifested  as  phase- wrapping  errors,  and  in  the 
latter,  e.g.  by  Marthi  V  Chengalur  (2014),  as  false  minima. 
We  have  demonstrated  that  in  either  case,  these  ambigui¬ 
ties  can  result  in  noticeable  distortion  of  the  reconstructed 
image. 

Using  the  Closest- Vector-Problem  formulation  of  the 
unwrapping  problem  due  to  Lannes  V  Anterrieu  (1999),  we 
have  developed  the  notion  of  a  wrap- invariant  pattern.  For 
wrap- invariant  patterns,  the  impact  of  the  2/r-periodicity 
can  be  completely  eliminated  (c.f.  Section  2)  using  well- 
known  algorithms  from  lattice  theory  and  the  Smith  Nor- 
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mal  Form,  and  reduced  to  a  mere  image  shift  (c.f.  Section 
4)  when  existing,  fast  approaches  are  used.  Phase-approach 
solutions  (Arnot  et  al.  1985),  (Greenaway  1990),  (Wieringa 
1992),  (Lannes  &  Anterrieu  1999))  are  commonly  used  to 
quickly  obtain  an  initial  point  to  aid  in  the  convergence  of 
the  Phasor  approach.  They  are  obtained  by  selecting  spe¬ 
cific  solutions  from  the  general  family  in  Equation  (24)  via 
enforcement  of  additional  constraints  on  the  solution.  Pha¬ 
sor  approaches  seek  those  complex  gains  and  object  visi¬ 
bilities  which  minimize  a  squared-residual  with  respect  to 
the  observed  complex  visibilities  (Wijnholds  &  Noorishad 
2012),  (Liu  et  al.  2010),  (Marthi  &  Chengalur  2014)  us¬ 
ing,  for  example,  gradient-descent  methods.  We  have  seen 
that  the  same  2/r- ambiguity  which  creates  a  family  of  Phase- 
approach  solutions  also  produces  a  corresponding  family  of 
global  minima  in  the  Phasor  approach.  For  either  case,  our 
results  show  that  for  wrap- invariant  patterns,  this  family 
represents  merely  shifted  versions  of  the  true  image.  We 
have  also  extended  this  analysis  to  show  that  wrap- invariant 
patterns  admit  reliable  imaging  using  standard,  and  general¬ 
ized,  phase  closures.  Conversely  we  show  that  patterns  which 
are  not  wrap- invariant  can  suffer  from  distortion  of  the  sort 
depicted  in  Figure  3. 

The  prognosis  for  mitigation  of  the  ambiguity  issues 
raised  in  this  paper  is  quite  positive.  Random  patterns  ap¬ 
pear  to  satisfy  the  wrap- invariance  condition  with  high  prob¬ 
ability.  Moreover,  failure  to  meet  this  condition  amounts  to 
the  existence  of  a  particular  kind  of  cycle  in  the  interfero¬ 
metric  graph  which  can  be  easily  isolated;  a  chief  contribu¬ 
tion  of  this  paper  is  a  simple  algorithm  for  identifying  such 
cycles  so  that  they  can  be  removed  by  the  array  designer. 
Finally,  we  have  shown  an  example  execution  of  the  algo¬ 
rithm  to  diagnose  a  member  of  a  popular  pattern  which  is 
not  wrap- invariant.  It  is  clear  that  with  careful  array  design, 
both  Phase-  and  Phasor-based  RSC  techniques  can  reliably 
produce  quality  image  reconstructions  free  from  discernible 
artifacts. 
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APPENDIX  A:  PROOF  OF  PROPOSITION  2.7 

In  this  Proposition,  we  propose  that  we  can  solve  for  the 
Fourier  phase  vector  from  the  closure  relations  in  two  sep¬ 
arate  integer-preserving  steps:  the  first  ( Step  A)  involving 
the  SNF  decomposition  of  Cm^c,  and  the  second  [Step  B ) 
involving  that  of  M.  To  show  the  validity  of  this  two-stage 
decomposition,  we  require  the  following  Lemma. 

Lemma  A.l:  ker( Cm_>c)  =  L 

Proof :  Clearly  every  vector  in  L  must  be  in  ker( Cm_>c), 
since  otherwise  there  would  exist  a  combination  of  {0}  which 
would  produce  a  non-zero  closure  phase.  This  is  impossible 
due  to  the  atmosphere-annihilating  property  of  the  closure 
phase.  We  must  also  show  that  a  basis  for  L  spans  the  kernel. 
Note  that  dim(im( Cm^c))  =  where  we  use  im(.)  to 

denote  the  range  of  a  matrix.  Noting  that  Cm_>c  has 

columns,  we  have  dim(ker( Cm^c))  =  Qf)  -  (A^_1)  =  N  - 
1  by  the  Rank-nullity  Theorem  from  linear  algebra.  Since 
dim(L)  =  N  -  1,  the  column  vectors  in  L  indeed  span  the 
kernel.  □ 

Suppose  we  have  the  SNF  decomposition  of  Cm_>c  = 
UCDCVC.  The  result  of  Step  A  is  then  a  particular  solution 
for  baseline  phases  given  by: 

Pp=VclV+c\S-c\y*cl+2xe*hci)  (Al) 

Given  Lemma  A.l,  Equation  (Al)  can  be  re-written  as: 


Substituting  Equations  (Bl)  and  (25)  into  Equation 
(26),  we  obtain: 

efr  =  V^E*  UjU^Vjfco  (B2) 

Noting  that  U &  is  orthogonal,  this  equation  can  be  sim¬ 
plified  to 

e<T  =  (Vo-  -  l\l)v£fc0  (B3) 

where  N  is  a  matrix  of  the  same  size  as  V^- ,  which  is  zero 
except  for  the  last  three  columns.  These  last  three  columns 
are  identical  to  those  of  V^-,  and  hence  by  Lemma  4.1  com¬ 
prise  an  orthogonal  basis  for  the  nullspace  of  M.  Noting  the 
orthogonality  of  V^-,  this  can  be  further  simplified  to: 

e<T  =  fco  -  NV^fco  (B4) 

To  proceed,  the  following  Definition  will  be  useful: 
Definition  B.l  (The  canonical  basis  for  the 
nullspace  of  M):  The  canonical  basis  {wi),i  e  1,2,3  for  the 
three-dimensional  nullspace  of  M  can  be  derived  trivially 
from  the  well-known  tilt-position  degeneracy  in  interferom¬ 
etry  described  in  Section  2.1.1  (Wieringa  1992).  Namely,  we 
can  define  the  basis  as  the  columns  of  a  (N  +  d)  X  3  matrix 
W^er(M)  as  follows: 


pp=^  +  vL+2nVc'D+cUc'elcl  (A2) 

where  J3q  is  a  phase  measurement  vector  corresponding  to 
the  optimum  least-squares  solution  for  Equation  (22),  and 
vp  is  a  vector  in  L.  Since  the  elementary  divisors  of  Cm_>c  are 
all  1  by  construction  (by  Lemma  2.6),  the  rightmost  term  in 
Equation  (A2)  will  be  0  mod  2 n  and  hence  the  first  step  is 
thus  integer-preserving.  We  then  consider  Step  B,  in  which 
we  solve: 


M 


e 

0 


= U mDmVm 


e 

0 


=h 


(A3) 


Since  \l  is  restricted  to  L,  its  impact  on  the  solution 
of  Equation  (A3)  is  restricted  to  a  perturbation  of  the  at¬ 
mosphere  phases  {0}  plus  a  linear  ramp  in  the  Fourier  phase 
(which  is  equivalent  to  an  image  shift).  If  all  elementary  di¬ 
visors  of  M  are  1,  application  of  the  SNF-based  generalized 
inverse  to  pp  will  preserve  the  integrality  of  the 

rightmost  error  term  in  Equation  (A2).  Hence  we  will  be 
able  to  recover  the  Fourier  phase  (up  to  an  image  shift)  in 
Step  B,  and  therefore  the  Proposition  holds.  □ 


APPENDIX  B:  PROOF  OF  LEMMA  4.2 

In  this  section,  we  prove  Lemma  4.2,  which  states  the  fol¬ 
lowing: 

Given  wrap -invariance,  the  (column)  vector  UxM^e^ 
has  integer  entries  below  row  2. 

Given  that  we  have  a  wrap- invariant  pattern,  we  know 
that  the  elementary  divisors  of  M  are  all  1.  Hence  there 
exists  an  integer  vector  Uq  such  that: 


e*h  =  Mfc0 


(Bl) 


W^r(M)  =  {wi\w2\w3}  = 


0  A  rx  A  ry 

A/Vxl  rx  ry 


(B5) 


where  is  a  length- A  vector  of  ones,  rx  and  ry  are  the 

x-  and  y-positional  coordinates  of  the  apertures  associated 
with  each  row,  respectively,  and  A rx  and  A ry  their  respec¬ 
tive  pairwise  differences.  □ 

Note  that  each  of  the  three  non-zero  column  vectors 
[vk},k  <=  1, 2, 3  in  N  can  be  expressed  as  a  linear  combinations 
of  the  elements  of  the  canonical  basis  {wi),i  e  1,2,3  defined 
above,  i.e. 


Vfc  =  a\w\  +  a2w2  +  <23^3  (B6) 

As  in  Section  4,  let  us  again  use  K  to  denote  the  set  of 
indices  in  eo-  associated  with  the  Fourier  phases  (as  opposed 
to  the  piston  phases),  and  their  corresponding  rows  in  N. 
Hence  we  have: 


Uxecr ,k  =  UXfco,Jf  -  UxIM^fco  (B7) 

Since  ux  is  an  integer  matrix,  the  first  term  is  clearly 
integral.  Let  us  then  examine  the  second  term,  and  in  par¬ 
ticular,  the  product  UxN^.  By  substitution  from  Equation 
(B6),  we  see  that  for  any  of  the  three  non-zero  columns  x 
of  we  have  (via  changing  the  basis  for  the  nullspace): 

Vxvk,K  =  aiVxwl,K  +  a2^Xw2,K  +  ai\Jxw3  ,K  (B8) 

where  {wj^}  are  the  vectors  comprising  the  upper  partition 
of  Equation  (B5).  The  first  term  in  Equation  (B8)  is  trivially 
0  since  that  w\  %  is  the  zero- vector.  Now  note  from  Theorem 
2.2  that  Ux  is  a  matrix  which  annihilates  all  the  spatial 
frequencies  in  the  matrix  X  below  row  2.  But  from  Definition 
B.l,  these  spatial  frequencies  are  in  fact  the  contents  of  the 
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two  columns  w^k  and  w$  k  (uP  to  a  uniform  scaling  factor). 
Therefore  the  column  vectors  in  Equation  (B8)  is  zero  below 
row  2.  This  means  in  turn  that  the  second  term  in  Equation 
(B7)  is  zero  below  row  2,  and  hence  that  is  integral 

below  row  2  (since  the  first  term  is  integral).  □ 

This  paper  has  been  typeset  from  a  TgX/ETEX  file  prepared  by 
the  author. 
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